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Abstract 

We consider a second order differential operator A(x) = — y\ dia i: j(x)dj + 

^2 =1 dj (bj (x) ■) +c(x) on R d , on a bounded domain D with Dirichlet boundary 
conditions on dD, under mild assumptions on the coefficients of the diffusion 
tensor ay. The object is to construct monotone numerical schemes to approx- 
imate the solution to the problem A(x)u(x) — fi(x), x g D, where \x is 
a positive Radon measure. We start by briefly mentioning questions of exis- 
tence and uniqueness, introducing function spaces needed to prove convergence 
results. Then, we define non-standard stencils on grid-knots that lead to ex- 
tended discretization schemes by matrices possesing compartmental structure. 
We proceed to discretization of elliptic operators, starting with constant diffu- 
sion tensor and ending with operators in divergence form. Finally, we discuss 
W^-convergence in detail, and mention convergence in C and L\ spaces. We 
conclude by a numerical example illustarting the schemes and convergence re- 
sults. 
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1 Introduction 

The object of present analysis arc numerical solutions of the elliptic boundary value 
problems in terms of monotone schemes. It is assumed that the elliptic differential 
operator has the divergence form with measurable coefficients satisfying the strict 
ellipticity condition. The basic idea of construction of monotone schemes is presented 
in [MWj without analysis of convergence of approximate solutions. Some elaborations 
of this basic idea can be found in [SMMM . LR2 and |LR3| where the convergence is 
considered in C- and Li-spaces. In these works the stencils of schemes are enclosed by 
rectangles with vertices at grid-knots. In the present work we extend the published 
results by constructing schemes with stencils stretching far from basic grid-rectangles 
and so being conceptually closer to the original idea of [MW . The schemes are not 
derived from finite difference operators approximating differential operators but rather 
from a general principle which ensures the convergence of approximate solutions. For 
the case of a classical elliptic problem this general principle is necessary and sufficient 
for the convergence in Holder spaces. 
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In Section [2] we describe boundary value problems for linear elliptic differential 
operators which are analyzed here from the numerical point of view. Problems for 
bounded domains and Dirichlct boundary conditions are our main interest. For the 
sake of completeness, numerical methods for problems on the whole real space are 
also studied. 

Monotone schemes for linear elliptic problems can be easily constructed by us- 
ing discretizations of the corresponding differential operators in terms of matrices of 
positive type. Therefore, in Section [3] we start our analysis by defining matrices of 
positive type and matrices with the compartmental structure. The latter ones are 
more pertinent to construction of discretizations of elliptic differential operators in 
divergence form. Monotone schemes give rise to grid-solutions. By embedding grid- 
functions into the space of hat functions we pass from grid-solutions to functions 
which are called approximate solutions. The remaining part of this section contains 
some technical results relating various norms of grid-solutions and the corresponding 
norms of approximate solutions. 

At the beginning of Section 2] we describe discretizations of the operator A = 
— ^ciijdidj in terms of matrices of positive type by using forward and backward 
finite difference operators [LRlj . Then we describe another method of discretization 
which is not based on any kind of finite difference operators. This method becomes a 
basis for discretizations of the elliptic operators in divergence form, A = — diO-ijdj. 
The resulting schemes are called extended schemes. The strong ellipticity of associ- 
ated discrctized bilinear forms [Yo is analyzed and proved for the class of extended 
schemes. 

The convergence of numerical solutions is analyzed in Section [5l First we con- 
sider the case of W\-cow&si%%ylcq and later the convergence in Holder spaces and 
Li-spaces. The convergence in H^-space is proved by using the standard finite ele- 
ment techniques. 

Some technical problems arising from a non-continuity of the functions ay, i =/= j, 
are discussed in Section [6l It is demonstrated by an example how the extended 
schemes can be applied straightforwardly to problems with general measurable func- 
tions aij. 

2 Preliminaries 

Elements of K d are denoted by x,w,r etc. The Euclidean norm in M. d is denoted by 
| ■ |. For a subset S C M. d the closure is denoted by S or cls(S) and boundary by dS 
or bnd(S). For an open set D the scalar product and norm of L 2 (D) are denoted by 
(•|-) and || • || 2, respectively. The norm of L p (D) is denoted by || • || p . The Sobolev 
Wp -spaces and -spaces are defined in a standard way [Mai [St] , In addition we 
need the Holder spaces C^{D),a e (0,1), defined as m 

[MalEi]- The convex set of 

positive Radon measures fi on 2} (D), J dfi — 1 (tight probability measures) is denoted 
by V(D). Then (v\fi) = j D v(x)fj,(dx) is well defined for v E W^(D). We say that 
a sequence of fi n S V(D) converges weakly to /i € V(D) if lim„(i;|^„ } = (v|/i) for 
each v € C(D). An open subset of R d with Lipshitz boundary [St] is called a domain 
with Lipshitz boundary. Only such domains are considered here. The following result 
is used. Let C^(D) be the linear space of continuous functions, having continuous 
partial derivatives and zero values on 3D. Then, for a domain D with Lipshitz 
boundary, the space C^(D) is dense in Wp(D), 1 < p < oo |Maj . 

Here we study elliptic operators on M. d in divergence form defined by 

d d 

A{x) = - J2 d i a ij (x)d j + ' ) + c(x), (1) 
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for which the coefficients are constrained as follows. The functions ay = dji, bi,i,j = 
1,2, ... ,d and c are measurable on M. d , c > and (x) converge to constant values as 
\x\ increases. Apart from this, the elliptic operator must be strictly elliptic, meaning 
that there are positive numbers M_, M, < M_ < M, such that the double inequality 



M\x\ 2 < 2_, aij(x)ziZj < M\x\ 2 , x G R d , (2) 



|2 



holds. The main part of A is denoted by A$. 

For sufficiently large and positive A there exist kernels (x,y) ^ KA(X,x,y) on 
D x D, such that (x,y) i— » \x — y\ d ~ 2 KA(X,x,y) are uniformly bounded on D x D 
and Holder continuous in x on the set D \ y. The integral operator with kernel 
Ka(\'>') is denoted by R(X,A). It maps measurable and bounded functions with 
compact supports into continuous functions on D and (XI + A)R(X, A) = I in C(D). 
Let £(£>) be the linear space of functions on D and = R(X, A)£(D). Then 

A(x) : D(A D ) i-> £(D) is denoted by A D . We say that is defined by A(x) on D 
with Dirichlet boundary conditions on dD. If £(-D) = L p (D), p < oo, the operator 
Ad is maximally closed from !D(^4d) c L p (D) into L p (D). 

Let us define a real bilinear form on W^D) x W p (D), l/p+ 1/q = 1, by: 



d 

E 



a „ 

a(v, u) — 2_] / a ij( x ) div(x) dju{x) dx 

&i(a;) div(x) u[x) dx + / c(a;) w(a;) u(a;) da;. 



(3) 



For a domain £> with Lipshitz boundary dD and for each pair v £ W^(D), u G 
Wp{D) n {Am G (Wl(D))^}, 1 <p <oo, there holds the Green's formula: 

a(v,u) = (u|Au). 

The Green's formula is also valid for each pair v € W^ C (D), u G Wl(D) n {Au £ 

no)}. 

The boundary value problem, to be studied in this work, is defined by 

A(x)u(x) = fi(x), x e D, . . 

u\dD = 0, W 

where D is a domain with Lipshitz boundary and the nonhomogeneous term n is either 
a positive Radon measure or /i 6 V7 2 -1 (-D). In the case of _D = R d we understand 
that the boundary condition in ([4|) is omitted. 

The variational formulation of Q for a solution u G W^(D), p — 1,2, has the 
following form: 

a(v,u) — (v\fi), for any v e Wq (D) . (5) 

Solutions of (HJ) and ([5]) are called strong and weak solutions, respectively. 
The following result is valid [LR2] . 

THEOREM 2.1 Let D be a bounded domain with Lipshitz boundary. For each p G 
[1, d/[d— 1)) there exists a unique weak solution u of (0) belonging to the class W p (D) 
and possessing the following properties: 

(i) There exists a positive number c depending on M_, M,p, D, such that the following 
inequality is valid: 

II" Hp, l < c - 



(it) If {(U ra : n G N} C V(D) converges weakly to a /i G V(D), then the corresponding 

!< 

strongly in L p (D) to u = A D fi 



sequence of weak solutions {u n : n G N} C W p (D), u n — Aj^[i n , converges 
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3 Grid-functions and W^-spaces 

Let the orthogonal coordinate system in R d be determined by unit vectors e^, and let 
us, for each n £ N, define a numerical grid G n on R d by vectors x — Y^,i=i h(n) fe^e;, 
where the grid-step h(n) is determined by h(n) = 2~ n . If not necessary the grid- 
steps h(n) are shortly denoted by h. For the sake of simple notation the grid step 
h is assumed to be the same for all coordinate directions. The obtained results of 
convergence are valid as well for grids with grid-steps depending on direction e,. 
To each v G G n there corresponds a grid-cube C(h,v) — ]J 1 [vj,Vj + h), where Vj 
are coordinates of v. Cubes C(h,v) define a decomposition of R d into disjoint sets. 
We say that the sets G n (D) — G n (~l D are discretizations of D. The grids G n are 
homogeneous with respect to translations in the direction of coordinate axes, i.e. for 
x £ G n and t = hpi e, ,pi £ Z we have x + 1 £ G n . There exist subsets of G n which 
are also homogeneous in the defined sense. Let r £ G n and r = (n, r%, . . . , r<j) £ N d 
be fixed. The set 

G„(r ,r) = {r + /t fc; r; e; : ki € Z} (6) 

is a subsets of G„. It is homogeneous with respect to translations for khiV^i, ki £ Z, 
grid-knots. There are vol(R) = Yit=i r i disjoint subgrids (O making a partition of G n . 
Each is denoted G n (R), where R stands shortly for 2d parameters rg e G n ,r £ N d . 
The index set of G n is denoted by /„. Similarly, I n (R) is the index set of G n (R). 

The discretization of a function u € C(R d ) on G„ is denoted by u„ and defined 
by values at grid-knots, (u n ) = u(x m ) where a; m = [mxh^m^h^ . . . ,rridh) 6 G„, 
and m = (mi, rri2, ■ ■ ■ , ma) is a multi-index. The function u„ is usually called a 
grid function. We denote the linear spaces of discretizations by l(G n ) or l(G n (D)). 
Elements of l(G n ) are also called columns. The corresponding L p -spaces are denoted 
by l p (G n ) or l p (G n (D)), and their norms by | • \ p . The duality pairing of v e l q {G n ) 
and u G l p (G n ) is denoted by ( v|u). The scalar product in l2(G n ) is denoted by ( - |- } 
and sometimes by (-|-). The norm of l p {G n {R)) is denoted by | • \r p . For p e [1, oo) 
we have 



|u| Bp = vol(R) "a-. 



l/p 



while Juj/foo = sup{|ufc| : k G I n (R)}. A matrix A n with indices corresponding to 
grid- knots of G n (or G n (D)) is said to be defined on G n (or G n (D)). We use the 
notation A n > if all entries of A n are non negative. 

The shift operator Z(x),x G M. , acting on functions / : R d <—* R, is defined 
by (Z{x)f\{x) = f(x + z). Similarly we define the discretized shift operator by 
(Z n (r, i)u n ) m = (u„) n , where n = m + rhei. 

Discretization of differential operators. A function u £ (7W(R ) has con- 
tinuous partial derivatives diU,i — 1,2,... .d. With respect to a grid step h, the 
partial derivatives are discretized by forward/backward finite difference operators in 
the usual way: 

Bi(t)u(x) = \{u(x + tei) - u(x)), d 

Bi(t)u(x) = f(u(x) - u(x-te t )). ' XizK ' ( '> 

Discretizations of the functions diU on G n , denoted by Ui(r)u n , V^(r)u„, are defined 
by: 

(Z7 i (r)u„) m = Bi(rh)u(x m ), (V i (r)u n ) m = 3i(rh) u(x m ). 

Then 

U t {r) = (rhy^Z^rJ) - I), 

V t (r) = {rh)- l (l-Z n {-r,i)) = U t (-r) = -C/,(r) T . 
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Therefore we have Ui(—r) = Ui(r) Z n (—r, i) — Z n (—r, i) Ui(r), and similarly for Vj(r). 

In accordance with the previous terminology, we say that di, Ylij di a ijdj e ^ c - are 
differential operators on R d or D. We say that their discretizations are defined on G n 
or G n (D). In particular, discretizations of the differential operator ([1]) are denoted by 
A n . Naturally, matrices A n are the main object in this work. We intend to analyze 
a class of discretizations A n of with compartmental structure. 

DEFINITION 3.1 (Compartmental structure) A matrix A = {dij}u is said 
to have the compartmental structure if it has positive diagonal entries, non-positive 
off-diagonal entries and positive or zero column sums: 

an > 0, dij < for i ^ j, ctj = fly > 0. 

A matrix with the compartmental structure with zero column sums is called conserva- 
tive. If A is compartmental then B = A T is called a matrix of positive type. A matrix 
of positive type B is called conservative if A — B T is (compartmental) conservative 
matrix. 

Hence, we consider a class of discretizations A n with the following properties: 

{A n )kk > 0, (A n ) k i < 0, I ^ fc, , s 

£ fc > 0. [ *> 

Additionally, there must exist a positive number a 2 such that 

cr 2 = sup h(n) 2 sup{(A n ) kk : kel n }. 

n 

Then Q n = a 2 1 — h 2 A n > and < cr 2 . This fact enables a decomposition of 

A n as A n — h~ 2 (o 2 I — Q n ). Hence, the resolvents 

h 2 1 °° 1 
R(X,A n ) = — - - A2/ff2 g (1 + A 2/ CT 2)fc Q*> 

are non- negative matrices on G n , \R(X, A n )\i < 1/A. 

Let /' C /. Then A' = {a^}//// is called a diagonal submatrix of A = {ay}//. A 
non- negative matrix A = {a y }// is called irreducible if for each finite index set I' C I 
the corresponding diagonal submatrix is irreducible. A non-negative matrix which is 
not irreducible is called reducible. Thus a non-negative matrix A is irreducible if for 
each finite subsets I' C I there exists m £ N such that the entries ptj of P = A m are 
positive for i,j G I'. 

LEMMA 3.1 Let A n possess either the compartmental structure or be of positive 
type. Then the matrix R(X, A n ), A > on G n is positive iff it is irreducible. 

Let the matrices A n on G n be discretizations of the differential operator |T]) and 
let there be defined linear systems: 

A n u n = jLt„, (9) 

where u„,/x„ 6 l{G n {D)). We say that the system ((9]) numerically approximates the 
boundary value problem |4]). The columns u„ are called grid-solutions. Obviously, fi n 
are discretizations of //. If the matrix A n is compartmental or of positive type we have 
A n = D n — B n , where the diagonal matrix D n has positive entries and non-negative 
matrix B n has zero diagonal entries. Finite difference equations (J9j> can be rewritten 
as u„ = D~ B n u n + D^ 1 [i n or componentwise as 

Uk = Y, { D n lB n) H n + dfefeMfe- 
l=ik 
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Hence, if A n is compartmental or of positive type, then the obtained numerical scheme 
is monotone |BSj . 

For a matrix A n on G n we define numerical neighborhoods 

M{x) = {yeG n : x = hk, y = hi, {A n ) kl ^ 0}. 

Obviously, numerical neighborhoods of a system matrix and stencils of the corre- 
sponding finite difference schemes are mutually related. 

Imbedding of grid-functions into K^-space 

Let us define the quadratic functional on l(G n ) by q(u) — |t/iU.|| and (?h(u) = 
2 R2 on l(G n (R)). It is understood that Qr = q for G n (R) — G n . 
There exist symmetric matrices Q n on G n such that <7_r(u„) = (u|Q„u). Discrete 
analogs of W^-spaces are spaces of those u„ 6 l(G n (R)) for which the norm | • \bs,i'' 

MW = |u|?q + te(u), (10) 

is finite. By convention | • (2,1 = | • \b.2,i for = 1. The subspace of grid- 
functions u € w\(G n {R)) for which u„ = Iq^d)^ is denoted by w\{G n (R,D)). 
Hence, w\(G n (R,D)) for r — 1 is denoted by w\(G n (D)). The restriction of qn 
on G n (R,D) is represented as gn(u) = ( u\Q n (D)u ), where Q n {D) is a symmetric 
matrix on G n (R, D). By using the fact that the negative Laplacean on C^ 2 '(D) and 
its discretizations — ^ Villi in l2{G n (D)) have positive minimal eigenvalues, we can 
derive the following result. 

LEMMA 3.2 Let D be bounded and Q n (D) be irreducible on G n (R,D). Then the 
norms U0\) and qn{-) x / 2 are equivalent in w\{G n (R,D)), 

to(-) 1/2 > P\-\m,x, 

where (3 is independent of n. 

Let us consider a norm | • | ^2.1 on l{G n {R)) defined by (fTU)) . Any such norm is a 
semi-norm on l(G n ). Our object of interest are quadratic functionals: 

| u |Lg,2,l = vol(R) ^2r I U Ir2,1 , ni , 

= NI + ^yEtilWu|| 2 < |u|| 2)1 . U£l{Gn) - (U) 

Then | • | a « g ,2.i is a norm on l(G n ). Unfortunately, it is not equivalent to | • 1 2. 1 
uniformly with respect to n. 

An element (column) u n G l(G n ) can be associated to a continuous function on 
M. d in various ways. Here is utilized a mapping l(G n ) 1— > C(R d ) which is defined in 
terms of hat functions. Let \ be the canonical hat function on R, centered at the 
origin and having the support [—1, 1]. Then z 1— » cf>(h,x,z) = \(h~ 1 {z — hx)) is the 
hat function on R, centered at x G M with support [x — h, x + h] . The functions 
z l— * 0fc(^) = Ili=i 4>{h,Xi, Zi),Xi = hki, define d-dimensional hat functions with 
supports S(l, x) — Ylii x i — h,Xi + h]. The functions </>&(•) G G„, span a linear space, 
denoted by E n (M. d ). Let u„ € l(G n ) have the entries u„fe = (u„)fc. Then the function 
u(n) = J2 ke % d u nk<f>k belongs to E n (R. d ) and defines imbedding of grid-functions 
into the space of continuous functions. We denote the corresponding mapping by 

: l(G n ) h-> E n (R d ). Obviously there exists <I>~ 1 : E n (R d ) h-» Z(G„) and the spaces 
l(G n ) and -EVi(R d ) are isomorphic with respect to the pair of mappings ^n,^" 1 . It 
is clear that E n (R d ) C E n+ i(R d ) and the space of functions U n E n (M. d ) is dense in 
L p (R d ),p G [1, 00), as well as in C{R d ). Let us mention that J2k <^ = 1 on K<i - 
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Now we consider another collection of basis functions. To each x — hk 6 G n (R) 
there is associated a d-dimensional hat function 



X 



— hki 
hri 



obviously, with the support S(r, x) = W^Xi — rih, Xi +rih]. They span a linear space 
denoted by E n (R,R d ). A gam we have Y,k = 1 on R d . The mappings 
cannot be applied to elements of l(G n (R)) and E n (R, R d ), respectively. Therefore 
we define restrictions <S> n (R) : l(G n (R)) -» E n (R,M. d ) and by the following- 

expression: 

u(n) = $ n (R)u n = { u n) k ipk- (12) 

h 

If we have to emphasize that u(n) is related to a particular set of parameters R, then 
we use an extended denotation u(R, n). 

The following results are necessary in our proof of the consistency in Section [SJ 

LEMMA 3.3 Let for each n e N there exist columns u„,v„ S l(G n (R)) such that 
ft- d |u„||j 2 < c, ^ d |v„||j 2 < c, where c > is independent of n. Let u(n) = $ n (i?)u„ 
and v(n) = $ n (i?)v„. If 

limh d sup \{Z n (k-l) -J)u„|| 2 = 0, 

" |fe-J| = |r| 



then 



Urn ^(v(n)\u(n)) - h d (v„ | u n ) R ) 



= 0. 



It is easy to check ||w(n)||! < ^ d | u ri| 2 - There holds a similar inequality for any 
pair diU, UiU n . First one verifies 

(div(n) | diu(n)) = (3i(nh)v(n)\ 3i(hn)(u(n)) = (Ui(n) v„ | U^n) u„), 

so that the following Theorem holds. 

THEOREM 3.1 Let sequences of functions v(n),u(n),n 6 N, be defined by U2\) . 
Then 



(ii) 



< h d min 



Eti (div(n)\d iU (n))\ < ^g R (v n ) 1 /2 gi j(u„)V2 

I^CrOvIa sup {|(ZKi)-7)U' i (r i )u| 2 : | 
|t/,(r,)u| 2 sup {|(Z(u;,i)-/)^i(r i )v| 2 : h 



Now we can get the following useful result involving the norm ||u(n) || 2j i and its 
averaged value: 



Hn)\\l M = \HR,n)\\l v 



(13) 



COROLLARY 3.1 Let u{n) — <I>„u„. There exists a 1 6 (0,1), independent of n, 
such that 

(l-a 2 )h d \n n \% 2A < \MR,n)\\l tl < h d \u n \% 2A , 

(l-<J 2 )h d \u n \ 2 avg ^ < \\u(n)\\l vgt2A < h d \n n \l vg ^. 
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An element u G W 2 1 (JS> d ) does not belong necessary to E n (R, M. d ). In order to 
approximate u with elements of E n (R,R d ) we define: 



H n ) = ll^fclli 1 OfcH^fc- 

k£l n (R) 



(14) 



The numbers ||V>fc|| 1 ~ 1 (V'fcl u ) are called Fourier coefficients of u. 

The basic result for our proof of convergence of approximate solutions is formulated 
by using the quantity T p (w, u) defined by: 



r p (w,u) 



The kernels 



E 



(Z(w) 
1 

||-0fc||l 



i>k(x) ipk(y) 



define an integral operator which is denoted by K n . Actually the kernels io n define a 
^-sequence of functions on 



d x M. d and K n converge strongly in L p -spaces to unity: 



COROLLARY 3.2 Let p G [l,oo]. Then 
ft) \\K n \\ p < 1. 

(ii) There is a positive number k(R) 7 independent of n, such that ||(7 — K n )u\\ p < 
k(R) sup{T p (w, u) : \u>i\ < hri\. 

(Hi) The operator K n 6 L(L2( 1 R d ), L,2(R. d )) has the spectrum equal [0, 1]. 

If necessary, we write K n (R) instead of K n to emphasize that the integral operator 
is defined with a basis ipk in E n (R,M. d ) for a particular value of parameter R. 

THEOREM 3.2 Let v,u G W%(M. d ) andu(n),v(n) be defined by [T$. Then 
(v(n)\u(n)) — (v\u) < c(R) min 

\\dv\\ 2 



\\uh sup lw]=h T 2 (w,v), 
II v h &vP\w\=h r 2 (w,u), 



(dv(n) \ dii(n)) — (dv\du) < c(R) min < 



where c(R) is n-independent. 



\\du\U 



\\du— 3 m 1 1 2 + sup T 2 (w,du) 

\W\=h 

\\dv — 3v\\2 + sup T2(w,dv) 

\W\ = h 



The procedure by which the norm (|TT|) is accompanied with the norms (ITU1) is 
thoroughly utilized throughout this work. Therefore we define the procedure gen- 
erally. Let a collection of functions u(R,n) defined by (fT2f or grid- functions u n (R) 
on G„(R) be denoted by £(R). For a functional F on €(R) or €(R) x €(R) the av- 
eraged value F avg is defined as the arithmetic mean of functionals over all or 
€(R) x €(R), respectively. For instance, the M^-norm on E n (R,M. d ) has its averaged 
version denoted by || ■ || Dg,2,i and explicitly expressed by (|13p . Another important 
example is a bilinear form a(-, ■) on E n (R,W d ) x E n (R, M d ). Its averaged version is 
defined by: 

a avg (v(n),u(n)) = — — - V a(v(R,n),ii(R,n)). 
vol(R) ^ 



4 Discretizations of elliptic operators 

Most of constructions and descriptions in this section are given for the differential 
operator Aq. At the end of section we briefly explain how to include discretizations 
for the lower order terms in (JTJ) - 
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4.1 Constant diffusion tensor 



Schemes and some of results regarding the convergence can be easily described for 
the case of classical second order elliptic operator Aa = —J2ij a ij^i9j, where the 
functions are Holder continuous on K d . At the beginning we analyze the case of 
constant diffusion tensor {atj} ff. 

A standard approach to a generation of discretizations of differential operators is 
based on utilization of finite difference operators approximating di, didj. So, by using 
the forward and backward difference operators Bj(pj/jj), 3 i{p%hi) we can define various 
discretizations 3uf(x) of didjf. In this way, the quadratic operator —aa{di) 2 — 
a,jj(dj) 2 is approximated by the standard central difference scheme aaUu + djjBjj 
with positive diagonal entries and non-positive off-diagonal entries. For i j finite 
differences are defined by: 

didjf (x) -> 3ijf(x) = 
f(x ± eihi ± ejhj) - f(x ± eihi) - f(x ± ejhj) + f(x), (15) 



hih 2 \ -f(x ± eihi =F ejhj) + f(x ± ah*) + f{x =p ejhj) - f(x). 

If ciij > 0, then aijdidj is approximated by the half sum of the first two possibilities, 
otherwise by the half sum of the second two possibilities. Hence, for the d-dimensional 

case 

Akkztr^ei = h' 2 ri ~ ^m^i rvT l^ m *l ' 

Akk±(r i e i +r j e j ) = ~ h 2 r . r . \ a ij\, a ij > 0> (^) 

^fcfc±(riei— rjej) h' 2 rir - \®"i3 I ' — ^' 



while the diagonal entries are calculated as the negative sum of off-diagonal ones. To 

Mj j 



a given diffusion tensor a — {dij}ff we associate an auxiliary tensor a defined by: 



HO — 



i ¥= J- 



By using Perron-Frobenius theorem one can easily prove the following assertion. 
There exist grid-steps h such that the matrix A n defined by (|16p has the compart- 
mental structure iff the auxiliary diffusion tensor a is positive definite. 

In case of d = 2 the tensors a, a are simultaneously positive definite or not. For 
higher dimensions a can be indefinite although a is positive definite. 

Apart from this standard approach to discretizations of differential operators there 
can be used constructions avoiding finite difference operators. One of possibilities is 
based on the following simple result: 

DEFINITION 4.1 A matrix A n with a finite band is called a discretization of the 
differential operator Aq = — a ijdidj if the equalities 

(Au)(k) = Y,( A n)klUi (17) 
I 

are valid for any polynomial x — > u(x) of the second degree. 

LEMMA 4.1 Let A n be discretizations of the differential operator Aq and u n = 
i?(A, A n )fi n such that \R(X, A n )\ oa is bounded uniformly with respect to n eFf. Then 
the sequence of functions il = {u(n) : n S N} C U n E n (D) converges in c^(m. d ) to 
the solution u — R(X,A )f of O). 

Proof: One can utilize arguments in standard proofs of the factor convergence 
of grid solutions as described in |RMj . In order to prove the convergence in C( a '(D) 
it is necessary to have the condition (fT7f as demonstrated in |LRlj . QED 
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A useful application of this lemma is the discretization defined by: 



Akk±ei 

A 



E 



A 



1 fcfc±( 
fefe±( 



h 2 r irj l a y I' 
' h' 2 r t r, \ a ij\> 



Tie.i—r;e.j 



a 



>o, 
<o, 



(18) 



where rj, rj G N. Of course, a convex combination of the system matrices (|18p is again 
a system matrix discretizing the classical differential operator Ao = — Y]^ aijdidj. 

In the case of two dimensions the obtained structures of system matrix can be 
classified into two groups by using the corresponding numerical neighborhoods. The 
grid G n — {h(kei + le 2 ) ■ k,l G Z} has the corresponding index set of indices 
k = (k,l). Possible numerical neighborhoods N{x) for the respective methods (Til 
(JTSJ), arc illustrated in Figure |4~T1 




(a): ai2 > 0, (b): a 12 < (c): a 12 > (d): a 12 < 

n = r 2 = i n = 3, r2 = l 

Figure 4.1 

In order to discretize (01 by monotone schemes we must have discretizations A n of 
fl} in terms of matrices of positive type. If we utilize discretizations (jT5)l the following 
inequalities must be satisfied: 



T- 



1 



, iOmi| < t iji, i € {1, 2, . . . , d}, 

flVm, tlTi 



and the following Lemma follows easily. 

LEMMA 4.2 There exist natural numbers fj such that the matrix A n defined by M8\ 
is of positive type iff the auxiliary diffusion tensor a is positive definite. 



As already mentioned, in higher dimensions a can be indefinite although a is 
positive definite. In such cases the simple schemes (fT8f cannot give us matrices of 
positive type. Rotations of coordinates are necessary to get system matrices of positive 
type. Naturally, so obtained schemes are also covered by Definition 14.11 

All the results of this subsection can be straightforwardly extended to cases with 
a general tensor- valued function x — * a(x). Numbers in Expressions (fTo| . (TTSJ) 
have to be replaced with the values dij{x) at the considered grid-knot x 6 G n . The 
condition of Lemma 14.21 must be satisfied uniformly with respect to x G G n . 



4.2 Operators in divergence form 

Discretizations of Aq(x) which are analyzed here are constructed by using several 
rules. Rules can be described vaguely as follows. The set M. d is discretized by 
G n {A) C G n , u by u„, diU by (Ui(pi)u n ), where pi are generally cc-dependent, 
are discretized by aij(xij(n,x)),x G G n (A), where xy(n, •) : K d i— > R d are certain 
functions, and integral is replaced by the corresponding Riemmanian sum of factors 
calculated at grid-knots. Thus the form is discretized by a sequence of forms 
v n ,u„ i — ► h(n) d ( v n \A n u n }. Without loss of generality, we may assume that entries 
of resulting matrices A n are defined by functions Xij(n,-) : R d ^ M. d and values 
aij(xij(n,z)),z G G n . 
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Standard schemes in two dimension 

We consider a subset of grid- knots G n (R) = r + hri e± + l hr 2 e 2 : k,l 6 zj C G % 
and the forms 



2 



a n (v,u) = (ai(/iri)u)(a;) a^a; + -(nei + r 2 e 2 )j (r3j(ftrj)«)(x). 

ij=l SBeG n (R) 

The non-trivial off-diagonal matrix entries of A n are defined as follows: 
(AOfefeinei = _ WT (n a ll + ~~ 77 l a 12 



fcfc±(riei — r 2 e 2 ) h 2 rtr 2 ' 12 I' 



and 



= a-ijfx + ~(±nei + r 2 e 2 )), = ay fx + -(±r x ei - r 2 e 2 ) 



If ai 2 < and brackets in (|19|) have positive values, then the scheme (|19|) defines a 
matrix A n with the compartmental structure; this scheme is called the first scheme. 
Numerical neighborhoods are defined by Af(~'(x) = x + AA '(0), where 

W (_) (0) = {0} U | ± hnei ± /ir 2 e 2 , ± (fin e x - ftr 2 e 2 )}. 

An illustration of a set AA~)(0) is given in Figure f4~TT b) . 

For the second scheme the grid G n {R) is defined as before while 

a n (v,u) = E?j=i 'ExeG n (R) (Bi((-) i ~ 1 hr i )v)(x) 

x dij(x + \{hriei - hr 2 e 2 )) (Bj((-) j ^ 1 hrj)u)(x). 

Off-diagonal matrix entries are: 

(AOfcfeinei = - hk^[ \.Tl a fi _ ~k I a i2 

(^")fcfc±r 2 e 2 = _ 7^7 (^ 2 " a 2 2 ± _ ^l"^!)' ( - 2U) 
(^™)fcfc±(r iei +r 2 e 2 ) = ~~ /i 2 nr 2 fl 12 • 



If ai 2 > and brackets in (|20[) have positive values the compartmental structure is 
again attained. The set G n (R) is unchanged while the numerical neighborhoods have 
different structure, Af^ + ~>(x) = x +Af( + \o), where now 

7V (+) (0) = {0} U |± hriet ± hr 2 e 2 , ±{hr l e 1 + hr 2 e 2 )\, 

as illustrated in Figure HTTT a). 

For each scheme there exists only one x^(n, ■) (instead of 3 mappings Xij(n, ■)) 
which is implicitly defined by expressions of &y • Thus, for the first and second 
schemes we have the following respective expressions: 



n,x,r) = x + h2 1 (riei+r 2 e 2 ), 
x^ + ' (n, x, r) = x + h2~ 1 (riei — r 2 e 2 ), 



(21) 



where a;*^) contain explicitly the parameter r. 
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. • • xi > In a, general two-dimensional case 

hn(iy» • cls(G n (l)) the function o 12 may vary and the 
compartmcntal structure of A n can 



be broken at some arid-knots. Be- 
1 , , sides, the function a 12 changes sign 

CtS[Lr n {2)) 



on R 2 and the compartmental struc- 

bnd(G n ( 1 ) ) ^ bnd(G n (2) ture cannot be realized by using only 

one scheme for all the grid-knots. In 

the case of continuous diffusion ten- 

. „ sor one can prove the following rc- 

Figure 4.2 * 

suit. There exists a partition M = 

U T U; such that the function a,\i does not change sign on Dj . To each 

there correspond grid-steps hrf^\l) such that the sets G^n\l) = n G„ are sub- 
grids and G n (R) = U T iG^ (I). To each grid-knot in int(G^\l)) we associate one of 
neighborhoods (x) and entries of A n defined by (fT9|) or (f20|) . Construction of en- 
tries of A n at grid-knots on the boundary of sets can bear a considerable amount 
of difhculties. For instance, one of simple cases is given by G n — G n (l) U G n {2) as il- 
lustrated in Figure l4T2l where we are faced with difficulties regarding the constructions 
of sets Af(x),x G {R 2 : x 2 = 0} and the corresponding entries of A n . 

Extended schemes in two dimensions 

The mappings x^ T \n, x, r) in R 2 are defined by (|2Tj) . where the superscripts (±) are 
related to the sign of a±2, i.e. to the cases a\ 2 > and a\ 2 < 0, respectively. The 
bilinear form on G n x G n 




aii(x { ^(n, x)) (Ei(h)v)(x)(3i(h)u)(x) 



i]\ x 



\n,x,r)) (Ui(rih)v)(x)(dj(rjh)u)(x) 

r. 



+ ^2 a ij^ x{ H n , x i r ))— (5i(h)v)(x)(Hi(h)u)(x) - (3 i (rih)v)(x)(Bi(r i h)u)(x) 

(22) 

can be written as 

a { n ] (v,u) = (v\A n u), 

where r £ N 2 , x^(n, x) = x^' (n, x, 1) by convention, v, u are columns with entries 
v(x),u(x), and A n is the resulting system matrix. By utilizing the standard varia- 
tional calculus we can derive entries of A n . In order to rewrite the entries of A n in a 
concise form we use the following abbreviations: 

at + (r) = a,ij(x + 2~ 1 (±r 1 he 1 + r 2 he 2 )), 
a if( r ) = a ij( x + 2^ 1 (±r 1 he 1 - r 2 he 2 )), 

andaff = aS 8 (l) by convention. For instance, afj + (r) — a.- L j (x^ (n, x, r)), a^ + (r) = 
ctij(x( + '(ri, x, r) — n/iei), etc. Now we have the following result: 

(^ n )fcfe±ei = 

=~fr (^-^I^WI), ( 23 ) 

(^ n )fcfc±(riei -r 2 e 2 
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Obviously, the constructed schemes are analogs of the schemes (fl8|) which are derived 
for the classical elliptic differential operator Aq = — J^.. ayc^c^. If a\i < and 

}£l{^-^tir)\} > o, «,A 7 ,fe {=,-}, 

the constructed matrix A n is compart mental. We call it the first extended scheme. 

The second extended scheme can be constructed analogously starting from another 
bilinear form. Let us define an (v, u) form a„ \v,u) by the following simple change. 
Each 3i(rih) in (j22|) has to be replaced with 3»((— )' /w$). The obtained form is 
denoted by a„ (w, u) and the corresponding system matrix A n has entries: 

( An )kk±e 1 = ~W i a ii ~ 77 I a i2 ( r )l) i 

(^We, = -jH^-i?l^(r)l), ^ 

(^™)fcfc±(r iei +r 2 e 2 ) = ~~ n^ 2 fl 12 ( r )' 

in a complete analogy with (I18j) . It is compartmental for ai2 > 0. 

It is important to point out that the schemes (|23|) . (|24|) are obtained by a guess 
from the schemes (fT5)l . i.e. they are a priori defined from some general principle 
which is not based either on difference equations or variational equalities. Then 
bilinear forms ([22)) are constructed so that variational equalities are equivalent with 
the schemes. This non-conventional approach in the construction of schemes has an 
apparent drawback. The forms (f2"2"| do not seam to be positive definite. Actually, 
their positive definiteness can be ensured in the | • | Q ug,2,i-norm, i.e. a weaker norm 
than | • |2,i-norm. This problem is analyzed in the next subsection. 

In the case of a general function ai2 on R 2 we cannot use only one of proposed two 
extended scheme. Here we give a construction of A n , which discretize the differential 
operator A , possessing the compartmental structure. Our intention is to construct 
a sequence = {ffly }if of tensor- valued functions converging *-weakly in L^R 2 ) 
to the original tensor- valued function (fT]), such that: 

1. Each a(') defines a form a^ r -*(w, u) and for any fixed pair v, u € W^H^ 2 ) 

a(v,u) = Una a^ r \v,u). 

r 

2. The forms (•, •) can be discretized by using described extended schemes. The 
resulting discretizations are denoted by h 2 an '(•, •). 

3. The matrices A^ defined by an\v,u) = {v\A n vl) have the compartmental 
structure. 

In the first step of construction we assume that the functions ay are uniformly 
continuous on R 2 . For each r G N we define the sets 

F r (+) = {xeR 2 : a l2 {x) > r" 1 }, 

F r (-) = {xeR 2 : a 12 (x) < -r" 1 }, r E N. 

F r (0) = M 2 \ F r (+) n F r (-), 

Due to the uniform continuity of functions ay there exists a positive e r such that for 
sufficiently large r and n: 

1. dist(F r {+), F r (-)) = e r > 0. 

2. The set F r (Q) has a finite number of connected components. 
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3. To each connected component of F r (0) there corresponds one connected set of 
F r (0)mnt(u xeGn C(h,x)). 

Now we define the sets 

G n ±] = lW r(± )nG„ C(h,x), 
Gl 0) = M 2 \ (GW UGH), 

and point out that for sufficiently large n the sets F r (0) and G^ ' have the same 
number of connected components. In addition, each connected component is equal to 
the union of interiors of certain cubes C(h,x) with vertices x s F r (0). In order to 
complete a partition of M 2 into subsets Gn we can extend the constructed sets G^ 
in the following way. Each non- trivial C(h, x) n G { ° ] can be attached to the nearest 
existing Gn . The first and second extended schemes have to be used for grid-knots 
in so completed Gn , respectively. The associated forms o^(u,u) are defined as in 
l|2"2"j) with the following difference. Instead of the sum over grid-knots x e G n , in the 
case of dn{v, u), we must have the sum over grid-knots in respectively. Then 

ah (v, u) — (v n \A n vin ) an d A n is compartmental. Of course, the parameter r varies 
on G n so appropriate notation is r(x). 

In the last step of the construction we assume that the functions a.y are mea- 
surable. We have to smooth the functions by a mollifier {\ n ■ n S N}, where 
Xn( x ) — n 2 Xi{nx) and xi is a non- negative function with support in the disc Bi(0). 
In this way we get a sequence of tensor- valued functions a(n) — a * Xn converging 
*-weakly in L^R 2 ) to the original tensor- valued function a. The previous steps 
of constructions must be applied to the differential operators with coefficients 
a VJ (n). 



4.3 Higher dimensions 

The index set of pairs 1(d) — {{ij} : i < j, i, j = 1, 2, . . . , d, i ^ j} has the cardinal 
number m(d) = d(d — l)/2. To each index {kl} € 1(d) we associate three coefficients, 

{kl} 1 {kl} 1 {kl} / ot -x 

a lk = — J a kk, a}i = — J Hi = a kh (25) 

and a bilinear form a^ kl ^(-, •), 

a^(v,u) = J a\f } (x)d i v(x)d j u(x)dx. (26) 

Apparently, for each v € W c j (lR d ), u 6 Wj^R ) with compact supports, the following 
equality is valid: 

a(v,u) = a {kl} (v,u). 
{kl}ei(d) 

To each of the forms a^ kl ^(-, •) we can associate a sequence of forms a,V^(-, •) and 
matrices A n k ^ constructed by using extended schemes. Then the matrix 

An = J2 A n l} > (27) 

{fcz}e/ 

is a discretization of Aq(x). In the construction of A n kl ^ , pairs r^ kl ^ (x) = (r\ kl ^ (x),rl kl ^)(x) e 
N 2 can be used, generally depending on x. 

In the remaining part of this subsection the objects of analysis are the compart- 
mental structure of (|27|) and strict ellipticity of forms v, u i— > (v|yl n u). 
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For x = hk G G n the neighborhoods Af(x) have structures which are mutually 
similar one to the other. If all <Xij(x) ^ then the minimal number of elements in 
jV(ce) is 1 + d + d 2 . In this case the set N(x) consists of its center, 2d-grid- knots 
on the ci-dimensional cross {±e^ : i = 1, 2, ...,d} and 2 grid-knots in each two- 
dimension plane. Generally, the number of grid-knots in two-dimensional plane may 
be larger than 2. Here we consider only the case of at most two grid-knots in two- 
dimensional planes. This demand has the following implication on the construction of 
discretizations An^ . Let the pairs e r , e s and e s , e t define two-dimensional planes and 
let An , An^ be the corresponding discretizations which are constructed by using 
parameters r^ rs \r^ st \ Then there must hold (r( rs ') s = (r^) s . In such case the 
off-diagonal entries of A n have the following general structure: 

(^")fcfc± e! 

i J ^ n )kk±m zj 

where an are certain convex combinations of values aa(x^ 1 ^), x\^ are certain points 
close to x which are defined by the extended schemes and rriij is either 7"je, + rjej 
or Ti&i — Tjej, depending on the sign of aij. 

Let us assume that the auxiliary tensor-valued function x i— ► a(x) = {dij(x)}ff 
fulfils the double inequality: 

M(a)\z\ 2 < ^ z i a ij (x)z j < 21(d) \z\ 2 , x<EM. d . (29) 

The following result is valid: 

PROPOSITION 4.1 Let be uniformly continuous on R d and let the double in- 
equality i29\) be valid. There exist a function x i— > r(x) G N d , to be used in the 
construction of matrices An kL \ such that has the compartmental structure. 

In order to prove this proposition one has to apply Perron-Frobenius theorem to 
the matrix a(x)~ 1 and construct the eigenvector z(x), corresponding to the minimal 
eigenvalue of a{x), with the property Zi(x) > m > 0, where m is ^-independent. 

If the coefficients a^ are not continuous, we have to use a smoothing procedure 
and apply the derived results to smeared out coefficients <Xij(n) = a^ * \n- 

Let us turn to the problem of strict ellipticity. A discrete form a n (-, •) on lo(G n (R)) x 
lo(G n (R)) is said to be strictly elliptic |Yo| if there exist two positive numbers 
M_(a n ),M(a n ) such that 

d d 

K(a n )J2 \U t (n)u\ 2 m < (u\A n u) R < M{a n )J2 \ U i(n)u\ 2 R2 . 

i=l i=l 

We shall say that a sequence of discrete forms a n (-, •) is strictly elliptic uniformly with 
respect to n G N if the positive numbers M(g„), M(a n ) are independent of n. 
The quantities an in (|28j) have the general structure 

oj n (a u ,x) = J2 s =i f s {x)a li (x li (n,x,s)), 

fs(x) > 0, Eti fs(x) = 1, 

where Xij (n, x, s) are defined by the construction of schemes. 

PROPOSITION 4.2 Let Aq = — ^ diOijdj be discretized by matrices A n with the 
following properties: 

a) Each A n has the compartmental structure whenever the double inequality &29\) 
holds. 



j_ 

'h 2 



—jp\ciij(Xij)\, 



(28) 
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b) There exist Pi,ti £ K such that 



{ A n) kk±r%ei = ~JfiP l w n (au,x) - ^ ti\a im (x im (n,x))\ 



at each grid-knot x = hk. 

c) If m is not in the direction of ±ej, i = 1,2, . . . , d, then {A n ) kk+m does not depend 
on any of an, % = 1, 2, . . . , d. 

Then the discrete forms v, u (■v\A n u) R are strictly elliptic on lo(G n (R)) x 
lo(G n (R)) uniformly with respect to n G N. 

Proof: Let us consider diffusion tensors a, b where a — {a^ }ff and b is defined 
by bij — aij — ndij. The corresponding auxiliary tensors are denoted by a, b as usual. 
Due to the fact that a, a are positive definite uniformly with respect to the points of 
R d there can be chosen n > so that b, b are also positive definite on R d . Let us 
define matrices H n on G n (R) by the following non-trivial entries: 

2k J 1 - , 1 K 
( Hn )kk = J^A^ ^2' ( Hn )kk±r ie , = ~^2> i=l)2,...,d. 

We have (v|.ff n u).R = k^^Li |^i( r i) u |?S2- The matrix B n — A n — H n has the com- 
partmental structure in accordance with assumptions of this proposition. A matrix 
B n with the compartmental structure is necessarily positive definite, i.e. there must 
hold (u\B n u)n > for any u e l (G n (R)). Therefore 

d 

(u\A n u) R = (u\H n u) R + (u\B n u) R > (u\H n u) R = K^\U l {r i )\i\ R2 , 

i=i 

proving the assertion. QED 
The first order differential operators of (TTJ) can be discretized by using upwinding 
method. In the case of continuous functions bj the resulting system matrices have the 
compartmental structure. If the coefficients bj are not continuous, the compartmental 
structure is obtained for sufficiently large n. 



5 Convergence of numerical solutions 

The difference schemes which are constructed in Section U enable us to generate grid- 
solutions u„ e L(G n (R)) and the corresponding imbedding u(n) — $„(i?)u„ into 
E n (R, M. d ) spaces. Basic schemes are derived from variational equalities 

( v„ | XI + A n u n ) = ( v„ | n n ), v„ e w\ (G„), (30) 

while extended schemes follow from a general principle which is not a priori related to 
any bilinear form and variational equalities. For a differential operators in divergence 
form Aq the derived discretizations A n have the compartmental structure, i.e. they 
satisfy (JTTJ) only for constant functions. Consequently, as expected, we cannot use the 
result on convergence of Section [3J In order to prove the convergence of approximate 
solutions we must replace conditions of Lemma 14.11 with some other property. It is 
natural to look for a formulation in terms of variational equalities pop and follow well 
known methods of finite elements. Unfortunately, we are faced with a difficulty. For 
extended schemes the bilinear form (f2"2"|) does not seam to be strictly elliptic. The 
strict ellipticity of forms is a basic supposition in a proof of -convergence. For this 
purpose we can use a general result of Proposition ^. 21 The underlining idea is clear 
from analyzing the 2-dimensional problem with a\i < on I 2 and a fixed r — (ri, r-i). 
When Proposition 14.21 is applied to the present case we get 
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COROLLARY 5.1 Let the auxiliary diffusion tensor a be strictly positive on R d . 
There exists a mapping x i— ► r(x) £ N d such that the matrices A n of have the 
compartmental structure, and the forms defined by a n (v,u) = (v|A„u) are strictly 
elliptic uniformly with respect to n £ N. 

For extended schemes Corollary 15.11 can be expressed simply as: 

2 2 

MJ2 [UiU n [l < a n {u,u) <MY1 \ U ^n\l 
i=i i=i 

with certain M, M. For basic schemes we must use averaged norms. In order to 
demonstrate that the averaged norms have no effect on the convergence of approxi- 
mate solutions we proceed further by using averaged norms. 



5.1 Consistency 

For constant functions ay we have an obvious equality: 

a(v,u) = lim h d a avg n (v(n),u(n)) 

n 

implied by (ii) of Theorem 13.11 and Corollary 13.11 For non-constat functions ay this 
property can be proved. Actually, we need a more general result. 

PROPOSITION 5.1 (Consistency) Let 2J = {v(n) : n £ N} and it = {u(n) : n £ 
N}, v(n),u(n) £ E n (R ), converge weakly and strongly in W 2 1 (R d ) to v,u £ W^R 4 ), 
respectively. Then 

a(v,u) = lim h d ( v„ | A n u„ ). 

n 

The present object of analysis is the functional F n (-, •) on E n (R, R d ) x E n (R, R d ), 
defined by 

F n (v,u) = h d vol(R) p(x(n,hk)) Bi(rih)v(hk) 3j(rjh)u(hk), 

kei n (R) 

where p is a bounded and uniformly continuous function on R d , x(n, •) is a trans- 
formation of R d such that \x(n, hk)) — hk\ < un" 1 , u> is fixed. Obviously, the form 
h d a n (-, •) of (|22|) is a finite sum of F avg:Tl with various p = aij. Therefore a proof of 
Proposition follows from the following result. 

LEMMA 5.1 Let the sequence = {v(n) : n £ N} C U n E n (R, R d ) converge W 2 l - 
weakly to v £ W 2 (R d ) and the sequence it = {u(n) : n £ N} C U n E n (R, R d ) converge 
W\ -strongly to u EW 2 1 {R d ). Then 

/dv du 
p(x) t—(x) - — (a;) dx. 
OXi OXj 

Proof: The integral on the right hand side is denoted by P(v,u). We have to 
analyze P(v(n),u(n)) as n — > oo. To simplify notation we assume i,j £ {1,2}. The 
indices k = (fci, k 2l ■ ■ ■ , k 8 ) are denoted shortly as k = (fci, k ) and k = (k\, k 2 , k ), 
k\,k 2 £ Z. The corresponding r±, r 2 are denoted by r,s, respectively. After inserting 
expressions for v{n),u{n) and carrying out a straightforward calculation, we get the 
following expression for the case of i = j = 1: 

P(v(n),u(n)) = ^2 dz(if} k '\p{z,-)'<Pl') v (.kk')iU(ki')i, 



k,k' ,1' £I„{R) 



7(k,r) 
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where v k i — Bi(rih)v(hk),Ukj — Bj(rjh)u(hk) and J(k,r) — [kh, kh + rh] . Analo- 
gously we get for i = 1, j = 2: 

P(v(n),u(n)) = dz 1 dz 2 (ip k "\p(z 1 , z 2 , 

k,l,h",l"ei n {R) Jj ( k > r )x J ( l > s ) 

X "(fc,i,fe'Oi U (W,H.J ^(^MO^) + W(t,I,fc")iM(Hr,I,(")j '/'Hr^O^fe) 
+ ^(fc,i+ s ,fc")i U (k,l,l")j ^k(zi)^l+s(z 2 ) + «(fc,J+s,h")i U(k+r,l,l")j ^k+r(zi)lpl +s (z 2 ) 

Now we approximate the obtained expressions by inserting specific values of p in 
the integrals over intervals J(k,r). The error arised from approximation can be 
estimated by a bound depending on the product of osc(p,n) = sup{|p(cc) — p(y)\ : 
x,y G K d , \xi — yi\ < 4tih} and the norms |{7j(rj)v|.B2, |tj( r i) u |-R2- For instance, 
for i = j = 1 the resulting expressing is 

Q n (v(n),u(n)) = h d vol(R) p{x(n, hk)) s k > r u^ V)j , 

k,k',l'eI n (R) 

where generally Ski — (V'fel^OIIV'fellr ■ The quantities sm have properties Ski > 
0, J2i s kl = 1- In addition 

|P(t!(n),u(n)) - Q„(u(rc.),u(n))| < osc(p, n) \v\ R2 \u\ R2 

< osc{p,n) (1 - a 2 )- 1 \\v(n)\\ 2 ,i ||tt(n)||a,i- 

Instead of P(v(n),u(n)) we consider in the following Q n (v(n),u(n)). The quantity 
Q n would be equal F n if s k i were absent and the double sum were replaced with the 
single sum over indices k. 

Q n (v(n),u(n)) - F n (v(n),u(n)) = h d vol{R) y~] p(x(n, hk)) vu Ski Sk'i' {uij-u kj ). 

kl 

To estimate the right hand side we need p = sup p. 

\Q n (v(n),u(ri)) - F n (v(n),u(n))\ < 
ph d / 2 \U i {r i )w n \ R2 h d ' 2 max{\(Z{ rj ,j) - l)Uj(rj)u n \ m : j = 1,2, . . . , d}. 

Due to the strong convergence of it, Lemma 13.31 and Corollary 13.11 we have for w = 
hr 3 ef. 

^K^r^i)-/)^^)^!^ < (l-o- 2 )-V 2 \\{Z(w)-l)u(n)\\ 2 
so that lim n F n (v(n), u(n)) = lim„ P(v(n), u(ri)) = P(v, u). QED 

5.2 W^-convergence 

It is described until now how the form X(v\u)+a(v, u) is discretized by forms h d X( v„ |u„ ) + 
( v n \A n u n }. In order to solve the discretized problem © or (|30[) . we have to describe 
how to discretize the linear function v — * ( v\fi) by h d ( v n \(i n ), where fi n € l(G n (R)). 
First we have to demonstrate the existence of /x n such that h d (v n \fji n ) — ► (v\fj,). 

LEMMA 5.2 Let \i be a continuous linear functional on W 2 (M. d ). There exists dis- 
cretizations fj, n (R) £ l(G n (R)) such that 



lim 

n 

for any R. 



h d (u n (R)\» n (R)) R - (u(R,n)\fi) 



\\u(R,n)\\, 



(31) 
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PROOF: It suffices to consider the case fi = df,fe L 2 (K d ), IHk-i = l/lk- Let 
us define 

(Mn) fc = -h- d (d^ k \f). 
By a straightforward calculation we get 

{u{R,n)\n) = ^ u k (rpk\df) = h d ^ u k Vk = h d (u\^ n ) R , 

k k 

proving Assertion. QED 
In other words the sequence of functions ji(n) = $„(i?)/i n (i?) 6 E n (R,M. d ) con- 
verges strongly in W^~ 1 (M d ) to fj,. This property implies the convergence of numbers 
( v(n)\n(n) } to the number (v\fi) for any W^-weakly convergent sequence of functions 
v(n). 

From the expression in proof of this Lemma we can get the estimate 
\{u n \n n ) a \ = hr d \(du(R,n)\f)\ < ft-« , |M| 2l _i||0u(J2,n)|| 2 < h~ d ' 2 \\nh,-i toK) 1/2 , 



where the last inequality follows from Theorem 13. II The obtained inequality is valid 
for any R, implying: 

KUftS.-l = Yl iMnlra,-! < ^ IH^-l- (32) 

Inequality p2|) and the variational equalities (|30j) imply the first result towards 
our proof of convergence of approximate solutions. If u„ solve (|30|) then 

Mn | avg,2,— 1 • (33) 

COROLLARY 5.2 Let u„ = R(\,A n )/j, n and u(R,n) = $ n (R)u n . Then for each 
R the sequence it = {u(R,n) : n G N} C U n E n (R, M. d ) converges weakly in W%(M. ) 
to some u G W 2 1 (R d ). 

Let it* be the solution of ([5]). Then the sequence of functions u*(n), defined 
by (|14p . strongly converges to it* in W%. In the remaining part of this analysis we 
have to demonstrate the expected property lim n u(R,n) = lim n u* (R, n) — u* for 
each R. We follow the well-known finite element technique. 

Mh d \u n -u* n [ 2 avgX1 < ^(u n -<|(A7 + A n )(u n -u;)) 
= h d {u n -u^\(XI + A n )u n ) 

- h d (u n -K\(XI + A n )K) 1 ' 

= h d (u n -u*|^ n ) - h d {u n -u^\(XI + A n )u*). 

By Lemma 15.21 the first term on the right hand side converges to (u — u*\ji). By the 
consistency property of Proposition 15.11 the second term converges to the same value. 



THEOREM 5.1 Let it be as in Corollary \5.SX Then the sequence it converges 
W 2 (M) -strongly to the unique solution u* of (5p. 

From this result, Lemma 15.21 and Lemma 13.21 we get another important result for 
A = 0. 

COROLLARY 5.3 Let D be a bounded domain with Lipsithz boundary and fi G 
W2 1 (D). Let A n (D) be the restriction to G n (D) of A n , // n on G n (D) satisfy [3l)) 
and u„ = An(D)^ 1 fi n . Then the sequence it converges strongly in (R. ) to the 
unique weak solution u of f5|). 
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5.3 Convergence in C- and Li-spaces 

Regularity properties of solutions of the converging sequences il of Corollarv l5.3l follow 
from the well-known results of DeGiorgi type. We can use criteria developed by jLU] 
as in our approach |LR3j in order to get the following result. If ||/i||oo = 1 then there 
exists a > such that il of Corollary 15.31 converges in the Holder spaces C^(D) to 
u of ([5]). Finally, if /i 6 Li(D), or p G V(D), numerical solutions must converge in 
the Banach space ti a \D), which is defined in [St]. This space is the completion of 
C^(D) in the norm: 

II "II = ||«||i + sup -i_||(Z(«;)-l) U ||i. 
\w\>o \ w \ 

Due to the fact that the spaces £^(-0,!^) are finite-dimensional, we can use the 
duality of finite-dimensional spaces li(G n (D)) and loc{G n (D)) and get a result on the 
convergence in Li(Z?)-space. In order to describe the convergence in Holder space and 
Li(D)-space we have to explain some details about discretizations. For the discussion 
of convergence in Holder space we assume p, = f € L QO (D), \\f\\oo = 1, while for the 
convergence in Li(D)-space we assume p G ^(D). Discretizations of / are denoted by 
f n and it is supposed that the sequence of functions f{n) — <I>f rl converges strongly to 
/ in L2(D). The grid solutions u„ = A^ 1 ^ define a sequence of approximate solutions 
it = {u(n) : n G N} G U n E n (D). Discretizations of p G P(D) are denoted by fi n and 
constructed so that p(n) G Li(D), lim„ ||/x(n)||i = 1, as well as lim„( g\p{n) } = (g\p) 
for each g G C(D). The corresponding grid-solutions v„ = A~ 1 /j, n define a sequence 
of approximate solutions 23 = {v(n) : n G N} G U n E n (D). 

COROLLARY 5.4 The sequence II has a subsequence converging to u in 

for some a > 0. The sequence 2J has a subsequence converging to u in (D) for 

some (3 > 0. 

6 Numerical realizations 

The convergence proofs of previous section were given for a simple case of only one 
r G N d realizing the compartmental structure of system matrices. In applications 
the coefficients vary and we can gain the compartmental structure by using 

several parameters r. The set M. d as well as D C M d must be decomposed into 
subsets Di which are defined during description of basic schemes. It is assumed 
that the compartmental structure of system matrices A n can be attained by a single 
parameter r(l) for each G n (l) — G n P\Di. Because of a numerous technical details, we 
cannot go into a thorough analysis of the construction of system matrices for variable 
r. Nevertheless, in the following example this approach is undertaken, and seems to 
be natural. 

Example 6.1 

In this example we consider the differential operator A = — &i a ijdj with the 

diffusion tensor defined as follows 

a(x) = pl Da {x), p 2 < a 2 , 

where a 2 is a positive number, p is a real number and Dq = (1/4, 3/4) 2 . 

Let D = (0, l) 2 G M 2 and dD be its boundary. The function x ^ u*(xi,x 2 ) = 
x\X2 is the unique solution to the following boundary value problem 

(Au)(x) = f(x) for x G D, 
u\dD = u*\dD, 



a 

a(x) 



a(x) 
1 
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where 



f{x) = 2pt Do (x) + (^(xi-^-Mfa- ~) + 6(x 2 - i) -35 fa- ~ 



The set R 2 is discretized by the grid G n of grid-knots x^i — hkei+hle2, k,l S Z, where 
/i is a grid-step. In order to get a discretization of D suitable for numerical handling 
we assume h = 1/N where N = AM. In this way we have discretizations G n (D) of the 
open square D = (0, l) 2 defined by grid-knots x^i = (xk, yi), 1 < k,l < N. The sets 
G n {2) = A) n G n and G n (l) = G n \ G n (0) define a partition of G n . The set G„(2) 
is closed, while G„(l) is open. Let G n {l,D) = G„(l) n D. Then G n (2),G n (l,D) is 
a partition of G n {D) to be used in constructions of numerical grids and approximate 
solutions. 

Figure 6.1 




To demonstrate the efficiency of the extended schemes we choose a 2 — 10, p = 2 
and the scheme parameters r\ — 1, T2 — 3 as illustrated in Figure [4~Tl part (c). These 
values of parameters ensure the compartmental structure of the system matrix. 

For a numerical illustration we discretize the domain D by 399 x 399 grid-knots 
with the grid-step h — 1/400, and define the grid-solution by the linear system ©. 
This system is then solved by iterations as follows. Let K n = diag(A„) and Q n = 
A n — diag(A„). Then Q n > and 

oo 

K 1 = E (K^Q n ) k K-\ (35) 

Let u„(r) be the approximation of u n after r iterations. By taking the stopping 
criteria to be |u„(r + l) — u„(r)|i < 10 -9 , we have found that the iteration terminates 
after r — 210 iterations. Then we compared the numerical approximations u„(210) 
and values of the solution in Loo(D) and Li(D)-norms. The L oc ,(_D)-norm of the 
difference u(hk) — (u ra (210))fc was 0.004, realized at the grid-knot with coordinates 
(0.75,0.75) which is one of the discontinuity points for a\2- The relative error in 
Li(Z?)-norm was estimated by 

_ |u* - U B (210)| 1 
E re l — 1UU j — -= , 

|u*|i 

where (u*), = u(hk). We obtained e re i — 0.0008, and the corresponding approximate 
solution u(n) as illustrated in the Figure ETT1 
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